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We study the self-averaging properties of the three dimensional site diluted Heisenberg model. 
The Harris criterion [l[ states that disorder is irrelevant since the specific heat critical exponent of 
the pure model is negative. According with some analytical approaches Q|, this implies that the 
susceptibility should be self- averaging at the critical temperature (R x = 0). We have checked this 
theoretical prediction for a large range of dilution (including strong dilution) at critically and we 
have found that the introduction of scaling corrections is crucial in order to obtain self-averageness 
in this model. Finally we have computed critical exponents and cumulants which compare very well 
with those of the pure model supporting the Universality predicted by the Harris criterion. 
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PACS numbers: 75.50.-Lk, 05.50.+q, 68.35.Rh, 75.40.Cx 



I. INTRODUCTION 

It is a well known fact that a possible way to obtain 
new Universality classes is to add disorder to pure sys- 
tems. The Harris criterion, Ref. [lj, states that if the 
specific heat diverges in the pure system, the disorder 
will change the critical behavior of the model, i.e a new 
Universality class will appear. Conversely, if the specific 
heat does not diverge in the pure system the critical ex- 
ponents of the disordered system will not change. 

The aim of this paper is to study the dependence of 
some observables with the disorder (self-averaging prop- 
erties) . The study of the self-averaging properties in dis- 
ordered systems has generated in the past years a large 
amount of both analytical, see Refs. {1, H, and nu- 
merical, see Refs. plla 0) works. We will focus in the 



averaging properties of the susceptibility both in vector 
and tensor channels. Results related to the universality 
of the critical exponents and cumulants are given in Ap- 
pendix A. Finally we give our conclusions in section V. 



II. ANALYTICAL PREDICTIONS 

The self-averaging (SA) of the susceptibility is defined 
in terms of: 



R\ 



(M 2 ) 2 - (M 2 ) 



(M 2 



(1) 



being M. the total magnetization. The susceptibility is 
self-averaging if R x — > as L — > oo. 



computation at criticality of the quantity R Y , which will T c ttj . , t , 

u„ La— a ,^ :„ .v,„„„i?' „„„„„„ In Ref - M, the following picture was found: 



be defined later and is a measure of the self-averageness 
of the susceptibility. 

We will study the three dimensional site diluted 
Heisenberg model with quenched disorder (in which, ac- 
cording to Ref. [H, the disorder is irrelevant) in order to 
test the analytical predictions. 

To obtain accurate measures of critical properties, in 
particular of R x , we will use Finite Size Scaling (FSS) 
techniques as the quotient method [1, [9[ that allows us to 
work in large lattices at the critical point and to perform 
infinite volume extrapolations. 

We will show results strongly supporting that this R x 
cumulant is zero at the critical point, but only taking into 
account the scaling corrections, against some theoretical 
predictions (Ref. 0]) and supporting others [3, H|. 

The structure of the paper is as follows. In the next 
section we summarize some analytical predictions con- 
cerning the self-averaging properties of diluted models 
at critically. In section III we define the model and the 
observables. In the first part of section IV we describe 
our simulation methods, in the second part we analyze 
deeply the correction to scaling exponent, u), while in the 
last part we give our numerical results concerning self- 



1. Outside the critical temperature: R x = 0. Can 
be found, based in Renormalization Group (RG) 
or using general statistical arguments, that R x oc 
(£,/L) d in a finite geometry, being L the system 
size and £ the correlation length, which is finite for 
T ± T c : then R x -> as L — >• oo. This is called 
Strong SA. 

2. At the critical temperature, a Renormalization 
Group analysis opens two possible scenarios: 

• Models in which according with the Harris cri- 
terion the disorder is relevant (a purc > 0): 
R x 7^ 0. The susceptibility at the criti- 
cal point is not self-averaging. In particu- 
lar, Ref. 2j shows that in these conditions 
R x is proportional to the fixed point value of 
the coupling which induces the disorder in the 
Hamiltonian, which controls the new Univer- 
sality class. This is called no SA. 

• Models in which according with the Harris cri- 
terion the disorder is not relevant (a purc < 0) : 
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0. The susceptibility at the critical 



point is self-averaging. In a finite geometry 
R x scales as L a l v — > 0, where a and v are the 
critical exponents of the pure system, which 
are the same in the disordered one. This is 
called Weak SA. 

The observable R x has been measured in other diluted 
models: for example in the four dimensional diluted Ising 
model, see Ref. |l0|. In this model a Mean Field com- 
putation and a numerical one found a non zero value for 
R x although the diluted model was shown to belong of 
the same universality class of the pure model, contradict- 
ing the conclusions of Ref. Q. One can claim that the 
logarithms which live in the upper critical dimension do 
the numerical analysis difficult. In particular was found 
analytically in Mean Field R x = 0.31024 and numeri- 
cally R x € [0.15,0.32]. Because the logarithms, it was 
impossible to do an infinite volume extrapolation for the 
numerical values of R x . Notice that in this model the 
only fixed point is the Gaussian one (all the values of the 
couplings are zero) and, following Ref. 0, R x should be 
zero. 

In addition a two-loops field theory calculation done in 
Ref. [U predicts a non zero value for R x for the diluted 
Heisenberg model (in which the disorder is irrelevant, 
"pure = —0.134, see Ref. [U). The two-loops field the- 
oretical prediction for a in the pure case was a pU rc > 0, 
so, apparently, this work is consistent with the findings of 
Ref. [H . The starting point in Ref. Q was the Mean Field 
computation done in Ref. [icj |. modifying it to take into 
account the vector degrees of freedom, introducing the 
fluctuations using the Brezin-Zinn- Justin (BZJ) method, 
Ref. 0. They found analytically R x = 0.022688 for the 
vector channel and universal (independent of the dilution 
for all p < 1). It is important to remark that in the BZJ 
method one fixes from the beginning the temperature of 
the system to the infinite volume critical one working in a 
finite geometry, so in order to compute R x in this scheme 
the following sequence of limits is used: 



R* = lim lim R X (L,T) 



T-*T C 



(2) 



where i?* is the infinite volume extrapolation at criti- 
cality of i? x (L,T), and T c is the infinite volume critical 
temperature of the system. The other possible limits se- 
quence that can be computed is: 



lim lim R^(L, T) , 

T^T c L~*oo X ' 



(3) 



which is zero even when the disorder is relevant since 
R x oc L~ d as T ^ T c . 

Hence, in order to test these discrepancies we have 
simulated numerically the site diluted three dimensional 
Heisenberg model computing R* in the vector and tensor 
channels. To perform this program, in particular in do- 
ing the infinite volume extrapolations of cumulants and 
exponents, is really important a proper use of the correc- 
tions to scaling. 



III. THE MODEL 

The Heisenberg site-diluted model in three dimensions 
is defined in terms of 0(3) spin variables placed in the 
nodes of a cubic three-dimensional lattice, with Hamilto- 



man: 



H 



-(3 ^ e i e j S i ■ S 3 > 

<i,j> 



(4) 



where the Si are three-dimensional vectors of unity mod- 
ulus, and the sum is extended only over nearest neigh- 
bors. The disorder is introduced by the random variables 
ei which take on value 1 with probability p and with 
probability 1 — p. An actual {ej} configuration, will be 
called a sample. 

In addition, as done in Ref. [8J, we define a tensorial 
channel associated with the vector S through the trace- 
less tensor: 



a/3 



S, 



1 



a,(3 = 1,2,3 



(5) 



In the following, we shall denote a thermal average 
with brackets, while the sample average will be over lined. 
The observables will be denoted with calligraphi c let ters, 
i.e. O, and with italics the double average, O = (O). We 
define the total nearest-neighbor energy as: 



£ — S ' Ci€jSi ■ Sj , 

(id) 



(6) 



and the normalized magnetization for both channels as 



Mf = yJ2^ S i 



?af3\ 



N 



(8) 



being V the volume (defined as L 3 , where L is the lattice 
size). Because of the finite probability to reach every 
minimal value for the free energy, the thermal average of 
Eqs. and (|5J), is zero in a finite lattice. Therefore, we 
have to define the order parameters as the 0(3) invariant 
scalars: 



M = {VM 2 } , M T = 
We also define both susceptibilities as 



(9) 



X = V(M 2 ) , XT = V(tvM 2 T ). (10) 
A very useful quantity is the Binder parameter, defined 



as: 



9i — 1 q 

3 (M 2 ) 



t , <(trM|) 2 ) 
.94 — 1 



(11) 
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Other kind of Binder parameter, meaningless for the 
pure system, can be defined as: 



v _ (M 2 ) 2 - (M 2 ) T {{tvM 2 T )) 2 - {tvM 2 T ) 



(M 2 ) 



<l2 



(trM|; 



(12) 

these are the quantities we are using to estimate the self- 
averaging properties of the susceptibility {R x ) in both 
channels. 

A very convenient definition of the correlation length 
in a finite lattice reads, see Ref. fl3i : 



x/f-i 

4sin 2 (vr/ J L) 



(13) 



where F is defined in terms of the Fourier transform of 
the magnetization: 



(14) 



V- 



F= —{\T{2tt/L, 0,0) | 2 + permutations) . (15) 
o 

The same definition is also valid in the tensorial case. 
This definition is very well behaved for the finite-size 
scaling (FSS) method we have employed, see Ref. 
Finally, we measure the specific heat as: 



about 10 5 of the first measures for L = 8 and growing 
this number with the lattice size. We have performed 
2 x 10 4 quenched disorder realizations in the diluted mod- 
els independently of the dilution and the lattice size and 
taking 100 measures per sample, according to Ref. [l(| 
who demonstrated that the best approach to minimize 
the statistical error is to simulate a great number of sam- 
ples with a few measures in each one. 

To measure the critical exponents, we use the so-called 
quotients method [1, [§] , which allows for a great statisti- 
cal accuracy. The starting point is the equation: 



Qc 



„xo/v 



0(L~ 



(17) 



With Qq being the quotient of the observable O mea- 
sured in a pair of lattices of sizes L and sL, at the temper- 
ature where — s, being u the eigenvalue correspond- 
ing to the first irrelevant operator on the language of the 
RG theory. We have fixed s = 2 in our case. Therefore, 
firstly we need to estimate by successive simulations the 
[3 point where: 



2L L 



(18) 



for each pair of lattices (L,2L), then we have used re- 
weighting techniques to fine-tune this condition. These 
re- weighting techniques are used to /3-extrapolate the ob- 
servables and calculate their /3-derivatives, always before 
the sample-average is performed. The equations used are, 
see Refs. 0,111: 



C = V- 1 {£ 2 ) - (£) 2 



(16) 



(0)(J3 + AP) = (Oe A ' 3£ )/(e A ^) 



(19) 



IV. SIMULATIONS 
A. Description 

The lattice sizes L we have studied are 8, 12, 16, 
24, 32, 48, 64, and, only in the pure model, L = 96. 

Between each measure of the observables described 
in Sec. firstly, we update the spin-variables using a 
Metropolis method over a 10% of the individuals spins, 
chosen at random, then we perform a number (grow- 
ing with L) of cluster updates using a Wolff method, 
see Ref. [lj], this is our Elementary Monte Carlo Step 
(EMCS). The number of clusters traced (or Wolff up- 
dates) between measures have been chosen to yield a 
good value of the self-correlation time, see Ref. [lj], in 
our case always 1 < t < 2 (r being the integrated auto- 
correlation time of the energy). 

In order to work in thermally-equilibrated systems we 
perform a great number of EMCS. We start the simula- 
tion always from random (hot) distributions of the spin- 
variables, although we have checked out that averages do 
not change if we begin from cold configurations (i.e. all 
spins being in the same direction). To concrete, we have 
taken 4 x 10 6 measures for the pure model discarding 



dp{0) = d p {0) = (0£-{0){£)). (20) 

These extrapolations are biased, for instance, the ex- 
pectation value of Eq. (|2"D|) , when the averages are calcu- 
lated with N m measures is: 



(' -£)%«». <*> 

hence, we have to correct this bias. We have followed the 
recipe given in Ref. [lOj. An example of the effect of this 
correction is found in Fig. Q] (for the greatest lattice size 
simulated andp = 0.9, the same applies to other dilutions 
and lattice sizes): a great bias affects the non corrected 
numerical data and it is clear the importance to take into 
account this effect. In addition, it is clear that the recipe 
of Ref. [l(| is working perfectly for N m = 100, which is 
the number of measures per sample we have performed in 
this work. Therefore, we are very confident that all the 
data presented in this paper, processing with the recipe 
of reference [lj| are not biased due to the reweighting 
technique. 

On the other hand, we have tried to use the solution 
obtained in Ref. [l5l |. where each sample is splitted in 
four parts, but the results were bad, this is due to the 
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FIG. 1: The g2 cumulant in both channels for L — 64, 
p — 0.9 with 1000 samples, /Simulation = 0.79112, reweighted 
at P = 0.79082 as a function of 1/N m , being N m the num- 
ber of measures in each sample. We report data with Nm = 
50, 100, 500, 1000, 5000, 10000. The data without the bias cor- 
rection proposed in Ref. [l(| are marked with triangles while 
the corrected ones are marked with circles. We also mark 
with the dotted lines the election used in this work (which 
corresponds to N m = 100). Notice the importance of the cor- 
rection of the bias if one performs reweighting with the data. 



small number of measures we take in each sample (10 2 ), 
which makes big differences between the averages in each 
quarter. 

To compute errors in the averages we have used Jack- 
Knife methods, see Ref. We have defined 20 Jack- 
Knife blocks for the pure model in a unique sample and 
one block for each sample in the diluted (p < 1) models. 

Calculated observables and critical exponents present 
sometimes, instead of a stable value, a monotonically de- 
creasing one. For n, such an evolution with growing L is 
found, but it is clearly weaker than for v. In these cases 
an infinite volume extrapolation is called for. If hyper- 
scaling holds, we expect finite-volume scaling corrections 
proportional to L~ w . This issue will be addressed in the 
next subsection. 



Exact results and RG calculations tell us that the dis- 
order, being irrelevant in this model, induce scaling cor- 
rections with an exponent ajv ~ —0.188 (in L) [17| . In 
addition to this new scaling correction one must have 
the one of the pure model, which is related to the cou- 
pling of the (0 2 ) 2 ) term in the Ginzburg-Landau theory. 
This exponent is assumed to be 0.8 []~8l.Tl9t] (for the pure 
model). Hence, the leading one is the exponent induced 
by the disorder. We will try to check this scenario by 
computing the "leading" correction to the scaling expo- 
nent from the numerical data. 

First of all, we have tried to estimate uj just by consid- 
ering it as another tunable parameter in Eq. (|22[) applied 
to some physical quantities. In these fits, as first approxi- 
mation, we disregarded the possible correlations between 
the data for different L values, the results can be seen 
in Table [T] If we perform a weighted average with this 
results we obtain uj = 1.07(9) for the pure model and 
u) = 0.92(9), oj = 0.81(7) and u = 0.88(4) for the diluted 
model with p = 0.9, 0.7 and 0.5 respectively, in very good 
agreement with the value from the references given be- 
fore. However, we think this method is not very confident 
because of the variability of the results from a quantity 
to another as shown in the Table |TJ 
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B. The scaling exponent uj 



As can be seen in Tables from IIVI to IXI| allocated in 
the appendix A, specially for the thermal exponents and 
the cumulants (g^ and 172), there are evident finite volume 
effects so we have to use the equation: 



V 



{L,2L) 



oc L 



(22) 



which is a consequence of scale-hypothesis, first derived 
in Ref. [IB]. Consequently, choosing a good value for uj 
is a primordial question. 



TABLE I: uj values from the L — > 00 extrapolations of some 
quantities. In the last row can be seen the weighted average 
of each column. We disregarded data with error bars bigger 
than the 100% of the values themselves. Those disregarded 
data are showed in the table as — . 



Another approach, following Ref. [20(, is to study the 
crossing points of scaling functions (as £/ ' L and (74) mea- 
sured in pairs of lattices with sizes L and 2L. The de- 
viation of these crossing point from the infinite volume 
critical coupling will behave as: 



A(3{L, sL) ee (3{L, sL) - /3 c (oo) oc 



1 



si -I 



-L~ 



(23) 
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With this method we need an additional estimate for 



the thermal exponent u 1 we have used, following llj, the 
value v = 0.7113(11) for the pure model (notice the really 
small error in is, so we will discard it in the following, see 
the comment at the end of this section), which is also a 
valid value for the diluted models, because of the validity 
of the Harris Criterion and as can be checked with the 
data in the Appendix A. Again we fixed s — 2. In this 
approach, we only use the crossing points in the vectorial 
channel because they are cleaner. 

Extrapolating this crossing points using Eq. (123|) . we 
can plot the minimum of the \ 2 of the fit as a function 
of uj obtaining the upper part of Fig. [2] and the whole 
Fig [3] To carry out this extrapolations we have to realize 
that the measures of the crossing points are correlated 
by pairs, so we have to use the x 2 definition that include 
the whole sclf-covariance matrix: 



for the diluted models with p = 0.9, 0.7 and 0.5 respec- 
tively, in a gree ment with the value obtained in the pure 
model pj. I 111 [HI, [l^], except in the p = 0.9 case, in 
which the value is far away two standard deviations from 
uj = 0.8 [H, Ojj]. One possibility is that we are comput- 
ing the leading correction to scaling exponent but with 
a large error. Another possibility is that in the p = 0.9 
model the coefficients of the leading correction (from the 
numerical point of view) to the scaling vanishes or are 
very small. This result and the change in the slope of 
the <?4 data for p < 1 respect to the p = 1 ones, as can 
be seen in Table UH are evidences of the possible im- 
proved action found for p = 0.9, see Ref. [15J, therefore 
the uj exponent we are measuring in this case could cor- 
respond to the third irrelevant operator, instead of the 
second one (remember that following RG the first one is 
a/v ~ -0.188). 



jV N 

X* = £ ~ fitHcov- 1 ),,™^ - fit), (24) 
1=1 m=l 

being N the number of crossing points, it is to say, the 
number of simulated L values minus two, xi is the ob- 
tained value for the observable x (£ y / L or gj), at the 
crossing point for Li and 2Li, and "fit" is the fitted value 
to the form of Eq. (f2"3"|) for L;. In addition: 



(cov) 



Lm 



\X m Xl) 



(x m )(xi) 



(25) 



can be also defined in terms of Jack-Knife blocks, see 
Rcf. ;ij], as: 



(cov; 



N h 



J-K 



,»• (26) 



where Nh is the number of Jack-Knife blocks, xj~ K are 
block variables, where the first subindex runs over L val- 
ues while the second one over J-K blocks, and (xi) is the 
average between all block variables given L = L\. 

Also, following Ref. [20fl . we can do a joined fit in u> of 
the crossing points of £ y /L and g\ by defining: 



^joined 



X\v /L 



■ X 



9l 



(27) 



using Eq. (|24[) to calculate each one of the latter terms 
and searching for the minimum of Xj 2 ; ne( j. We can ob- 
tain the error in uj searching for the point uj\ in which 
xfoincdK) = xfoincdKun) + 1 , so the error is Auj = 
l^min — The results for this joined fits are shown in 
the upper part of Fig. [5] and in the whole Fig. [3] We find 
with this method the values: 

uj = 0.96(15), 

for the pure model and: 

uj = 2.29(70) , 0.84(17) , 0.64(13), 
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FIG. 2: Upper part: \ 2 as a function of w deduced from the 
fits to L — * co, Eq. (|23[) . for the crossing point of £/L and qa 
for the (L, 2L) pair for the pure model. Also is shown the 
joined x 2 - Lower part: extrapolated /3 c (oo) as a to function, 
the point where both observables give the same extrapolated 
value is marked with the dotted line. 



In addition, as Ref. [20( also did, we could estimate the 
correct value for uj as the one producing the same /? c (oo) 
value for £/L and (74, as can be seen in the lower part of 
Fig. [21 marked with the dotted line in uj = 0.88. This 
approach only works for the pure model in which such a 
point is found, with another p value the /3 c (oo) estimates 
from £/L and 54 do not cross each other. 

In conclusion, we have shown that our data (both for 
the pure and diluted model) are fully compatible with the 
value uj = 0.80(1), obtained previously both numerically 
and analytically for the pure model [21| . In addition, 
since the error bars in uj are really small (one per cent 
of error) we have discarded the uncertainty in uj in the 
analysis presented in this paper, the error bars in the 
extrapolated quantities are much bigger than the uncer- 
tainty caused by the error bars in w: so we have fixed 
uj = 0.80. The extrapolations (Tables III to XIII) and 
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FIG. 3: Upper part: x as a function of for the diluted 
(p = 0.7) model. Also is shown the joined x 2 - Lower part: \ 2 
as a function of ui for the diluted (p — 0.5) model. 



Figures (3 to 8) shown in the rest of the paper are ob- 
tained using this value of to. 

Finally, it is interesting to note that we have seen in 
the analysis presented in this subsection no traces of the 
leading correction to the scaling exponent even for the 
strongest dilution we have simulated, which should be 
a/v ~ —0.188. One can explain this fact assuming that 
the amplitudes of this scaling correction exponent are 
really small, so we are seeing the next to the leading 
scaling correction. 
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0.1130(11) 


0.2061(24) 


0.6006(4) 


0.4999(8) 




12 


0.0834(8) 


0.1600(18) 


0.6072(3) 


0.5047(6) 




16 


0.0702(7) 


0.1395(16) 


0.6107(3) 


0.5070(6) 




24 


0.0553(6) 


0.1138(13) 


0.6138(2) 


0.5085(5) 




32 


0.0474(5) 


0.0980(11) 


0.6151(2) 


0.5095(4) 



TABLE II: Cumulants for the 0(3) model. In the first col- 
umn is represented the spin density p. All the cumulants are 
calculated in the crossing points of £/L for L and 2L. The 
averages have been computing using 10 4 samples (except in 
the p = 1) case. 



C. Numerical Results about Self- Averaging 

Once that we have checked that the value u> = 0.80 
describe the corrections to the scaling in the pure and 
diluted model we can try to extrapolate the values of the 
two g 2 to infinite volume. 

Numerical results for g 2 and (74 in both channels are 
shown in Table |TT] both for pure (only g±) and diluted 
models. 

First of all, we will try to check the non zero 172 scenario 
with the correction to the scaling exponent fixed to that 
obtained in the previous section. 

We have found that it is possible to extrapolate us- 
ing the form of Eq. (|22[) (performing a joint fit [22j) the 
values of g 2 to a value (depending only on the channel) , 
which is independent to the dilution, and near the ana- 
lytical prediction of reference [3J|. However, simulations 
(with a small number of samples), at dilutions p — 0.95 
and p = 0.97 do not follow the scaling found for p < 0.90 
(see Table HIT]) . Hence, our numerical data do not sup- 
port the scenario g 2 ^ 0, see Figs. 0] and Notice, see 
also Table IIII1 that all the values for these two lowest di- 
lutions are smaller than the extrapolated point and they 
are decreasing (for both channels and taking into account 
the error bars). 



Secondly, we will check the g 2 — scenario. To do 
this, we extrapolate 172 using the form proposed in Ref. 
(52 ^ L a l v ) but also including the term L~ w : 

32 = aL a ' v + bLT u . (28) 

We obtain the fits shown in Figs. [6] and [7] for both chan- 
nels. The x 2 of these all fits are really good, hence, we 
have obtain strong evidence supporting this g 2 — sce- 
nario. Notice that the introduction of the scaling correc- 
tion has had paramount importance in order to obtain 
very good x 2 m an the fits. The numerical data, for 
the lattice size simulated, does not follow the one term 
dependence g 2 cx L a ' v . 

V. CONCLUSIONS 

We have studied the critical properties of the Heisen- 
berg diluted model for different values of the dilution 
using the quotient method. Our main aim was to check 
the self-averaging properties of the susceptibility. 

We have studied in a great detail the corrections to 
the scaling in the diluted Heisenberg model. We have 
obtained that the numerical data follow the next to the 
leading correction to the scaling exponent instead the 



7 



V 


L 


gX 


T 
92 


gX 


T 


0.97 


8 


0.0108(6) 


0.0181(13) 


0.6201(4) 


0.5187(10) 




12 


0.0102(6) 


0.0189(14) 


0.6195(4) 


0.5164(10) 




16 


0.0084(6) 


0.0158(12) 


0.6201(4) 


0.5159(10) 




24 


0.0072(5) 


0.0146(11) 


0.6206(4) 


0.5162(9) 




32 


0.0074(5) 


0.0152(12) 


0.6206(4) 


0.5152(10) 


0.95 


8 


0.0179(10) 


0.0290(18) 


0.6180(5) 


0.5158(11) 




12 


0.0167(9) 


0.0329(20) 


0.6182(5) 


0.5116(12) 




16 


0.0150(9) 


0.0286(18) 


0.6181(5) 


0.5129(11) 




24 


0.0117(7) 


0.0228(14 


0.6186(4) 


0.5135(11) 




32 


0.0118(7) 


0.0251(17) 


0.6193(4) 


0.5140(10) 



0.25 



0.15 



TABLE III: Cumulants for the 0(3) model with high p values. 
All the cumulants are computing averaging only 10 3 samples. 




0.05 



FIG. 5: Joined extrapolation to L — *■ oo for the (?2 cumulant 
of the tensorial susceptibility, to the form of Eq. (|22[) . 




FIG. 4: Joined extrapolation to L — > oo for the gi cumulant 
of the vectorial susceptibility. Extrapolations are carried out 
by choosing a common value for the first term of Eq. (I22p 
for all dilutions and by minimizing the joined \ 2 ■ We have 
disregarded the data with L — 8 to obtain a good value of the 

x 2 . 



0.08 
1/L 

FIG. 6: Extrapolation to L — » oo for the gi cumulant of the 
vectorial susceptibility. The fitting function is in this case of 
the form g 2 = aL a/u + bL~ u '. 



leading one. We will show in the appendix all the critical 
exponents and cumulants using this next to the leading 
exponent, also we report that the result of this analy- 
sis is fully compatible with the RG predictions and the 
Harris Criterion: our exponents and cumulants are com- 
patibles with that of the pure model and independent of 
the dilution with a high degree of precision. 

In addition, we have shown that we obtain non uni- 
versal quantities if we assume a/v as the main scaling 
correction even if we add the ui correction to the scal- 
ing exponent (see the appendix) , using two correction to 
scaling exponents in the analysis. 

Finally, we have shown strong evidence for a zero gi 
cumulant, both in vector and tensor channels, in the ther- 
modynamic limit at critically, contrasting with some an- 
alytical predictions [3(], and in agreement with the ones 



obtained in Ref. 0. The introduction of scaling correc- 
tions in the analysis has become crucial to obtain the 
#2 = scenario. In addition, simulations of samples with 
lower dilution have helped to discard the g^ ^ scenario. 
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FIG. 7: Extrapolation to L — > oo for the g2 cumulant of the 
tensorial susceptibility to the form 172 = aL a ' u + 6L~". 
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APPENDIX A: CRITICAL EXPONENTS AND 
CUMULANTS 

In this appendix we will check the consistency of the u> 
exponent computed in the text by means of the compu- 
tation of critical exponents and cumulants. In addition, 
we will check if these sets of exponents are Universal or 
not by comparing different dilutions with the pure model. 
We will use in this analysis the data from p = 0.9, 0.7 and 
0.5 (we have simulated in these values of the dilution 10 4 
samples) . 

Eq. (fT71) applied to the operators dpi;, dpg/t, M and 
X, yields respectively the critical exponents 1 + X/v, 
(d — 2 + rj)/2 and 2 — 77. The numerical results are shown 
in Tables|TV]and|V]forthe pure model, Tables EJ and EED 
for p = 0.9, Tables Enll and US for p = 0.7 and Tables IS 
and IXII for p — 0.5. We have also carried out a joined 
extrapolation for every p values by fixing the same value 
of the extrapolated exponent (first term in Eq. I|22p) for 
every p value and then minimizing the joined \ 2 . Some of 
these fits are shown in figures [8] to QT] and the compared 
results can be seen in Tables IXIII and IXIII1 

The joined extrapolation of the Binder cumulant 34 is 
shown in Table IXIVI and the agreement of our results 
with the ones obtained in Refs. (numerical for the 
pure model) and [l| (analytical) is really very good. We 
obtain also complete agreement with previous numeri- 
cal estimates of the pure model critical exponents, see 
Ref. [HI- 

We obtain non universal critical exponents and cumu- 



lants if instead u> = 0.8 we use u = —a/v as the correc- 
tion to scaling exponent. In addition, the dilution depen- 
dent exponents and cumulants are clearly different of the 
pure ones. Furthermore, this scenario does not change if 
we fit the data using both oj = —ajv and u> = 0.8. 





n 


Vt 


L 


X 


M 


XT 


Mr 


8 


0.0301(7) 


0.0319(8) 


1.4301(12) 


1.4343(13) 


12 


0.0339(7) 


0.0353(8) 


1.4324(11) 


1.4352(12) 


16 


0.0348(7) 


0.0358(8) 


1.4310(11) 


1.4335(12) 


24 


0.0361(6) 


0.0367(7) 


1.4293(9) 


1.4307(10) 


32 


0.0369(7) 


0.0374(7) 


1.4289(11) 


1.4300(12) 


48 


0.0373(6) 


0.0378(7) 


1.4271(9) 


1.4280(10) 


L — > 00 


0.0391(9) 


0.0390(10) 


1.4250(13) 


1.4249(15) 


X 2 /d.o.f 


0.138/3 


0.354/3 


1.047/3 


1.952/3 


prob 


0.987 


0.950 


0.790 


0.582 



TABLE IV: Magnetic exponents for the pure 0(3) model. 
The last three rows correspond to the L — > 00 extrapolation 
(disregarding data with L = 8) 





V 


L 






dpgj 




8 


0.7016(30) 


0.7217(13) 


0.6846(41) 


0.7306(14) 


12 


0.7033(32) 


0.7162(14) 


0.6931(49) 


0.7188(13) 


16 


0.7028(35) 


0.7123(16) 


0.6830(56) 


0.7118(17) 


24 


0.7061(37) 


0.7123(17) 


0.6908(47) 


0.7112(18) 


32 


0.7081(35) 


0.7121(19) 


0.7022(61) 


0.7116(23) 


48 


0.7101(41) 


0.7118(19) 


0.7125(61) 


0.7085(21) 


L — > 00 


0.7109(38) 


0.7071(19) 


0.7082(51) 


0.7071(35) 


X 2 /d.o.f 


0.667/4 


4.104/4 


7.039/4 


0.565/2 


prob 


0.954 


0.392 


0.134 


0.754 



TABLE V: Thermal critical exponents for the pure 0(3) 
model. In the last column we have disregarded data with 
L < 16. 
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»7 


Vt 


L 


X 


M 


XT 


M T 


8 


0.0346(26) 


0.0345(28) 


1.4154(36) 


1.4176(37) 


12 


0.0360(24) 


0.0360(26) 


1.4195(34) 


1.4207(36) 


16 


0.0371(23) 


0.0374(25) 


1.4207(34) 


1.4218(35) 


24 


0.0373(22) 


0.0375(24) 


1.4204(32) 


1.4221(34) 


32 


0.0383(21) 


0.0383(23) 


1.4219(31) 


1.4227(33) 


L — > oo 


0.0397(29) 


0.0399(31) 


1.4245(41) 


1.4252(43) 


X 2 /d-o.f 


0.292/3 


0.124/3 


0.544/3 


0.137/3 


prob 


0.962 


0.989 


0.909 


0.987 





V 


L 


dpgY 




dfjgl 




8 


0.7888(69) 


0.7881(31) 


0.8256(143) 


0.8422(42) 


12 


0.7810(74) 


0.7806(33) 


0.8078(140) 


0.8067(41) 


16 


0.7633(70) 


0.7760(35) 


0.7739(131) 


0.7897(43) 


24 


0.7491(66) 


0.7628(37) 


0.7719(146) 


0.7792(47) 


32 


0.7400(67) 


0.7521(42) 


0.7656(178) 


0.7627(56) 


L — > oo 


0.7206(88) 


0.723(10) 


0.729(19) 


0.7255(61) 


X 2 /d.o.f 


2.313/3 


0.281/1 


1.314/3 


1.965/3 


prob 


0.510 


0.596 


0.726 


0.580 



TABLE VI: Magnetic exponents for the diluted 0(3) model 
with p=0.9. Extrapolations are carried out without disregard- 
ing data. 



TABLE IX: Thermal exponents for the diluted 0(3) model 
with p=0.7. In the second column the fit is obtained disre- 
garding data with L < 16. 





V 


L 


dpgX 








8 


0.7319(49) 


0.7443(24) 


0.7128(83) 


0.7709(29) 


12 


0.7381(53) 


0.7411(25) 


0.7267(86) 


0.7514(29) 


16 


0.7430(55) 


0.7381(26) 


0.7536(99) 


0.7426(31) 


24 


0.7384(57) 


0.7368(28) 


0.7337(95) 


0.7395(32) 


32 


0.7398(54) 


0.7365(29) 


0.7241(97) 


0.7345(33) 


L — > oo 


0.734(15) 


0.7318(33) 


0.728(17) 


0.7152(39) 


X 2 /d.o.f 


0.134/1 


0.168/3 


5.468/2 


3.156/3 


prob 


0.714 


0.983 


0.065 


0.368 



TABLE VII: Thermal exponents for the diluted 0(3) model 
with p=0.9. In the first and third columns we obtain bad 
results because not monotonally-decreasing series. 





V 




L 


X 


M 


XT 


M T 


8 


0.0505(45) 


0.0461(48) 


1.3435(61) 


1.3431(62) 


12 


0.0448(39) 


0.0439(42) 


1.3684(54) 


1.3702(56) 


16 


0.0421(36) 


0.0417(39) 


1.3877(51) 


1.3896(52) 


24 


0.0396(32) 


0.0406(35) 


1.4033(46) 


1.4053(48) 


32 


0.0399(30) 


0.0414(32) 


1.4126(43) 


1.4152(45) 


L — > oo 


0.0346(60) 


0.0378(46) 


1.446(12) 


1.449(12) 


xVd-o.f 


2.225/2 


2.191/3 


0.119/1 


0.327/1 


prob 


0.329 


0.534 


0.730 


0.568 



TABLE X: Magnetic exponents for the diluted 0(3) model 
with p—0.5. In the third and fourth columns we have only 
used data with L > 12. 





n 


TjT 


L 


X 


M 


XT 


M T 


8 


0.0436(38) 


0.0412(41) 


1.3882(52) 


1.3879(53) 


12 


0.0411(34) 


0.0401(36) 


1.4005(48) 


1.4007(49) 


16 


0.0392(31) 


0.0392(34) 


1.4061(45) 


1.4073(46) 


24 


0.0383(29) 


0.0386(31) 


1.4131(41) 


1.4136(43) 


32 


0.0382(27) 


0.0389(29) 


1.4142(40) 


1.4149(41) 


L — > oo 


0.0343(57) 


0.0370(58) 


1.4299(72) 


1.4318(76) 


X 2 /d.o.f 


0.232/3 


0.059/3 


0.472/3 


0.567/3 


prob 


0.972 


0.996 


0.925 


0.904 



TABLE VIII: Magnetic exponents for the diluted 0(3) model 
with p=0.7. Extrapolations are carried out without disregard- 
ing data. 





V 


L 


dpgY 




dfjgl 




8 


0.8102(91) 


0.8357(46) 


0.9180(241) 


0.9540(72) 


12 


0.8042(90) 


0.8322(50) 


0.8880(248) 


0.8866(71) 


16 


0.7764(89) 


0.7862(48) 


0.8449(242) 


0.8136(64) 


24 


0.7702(93) 


0.7778(52) 


0.8311(234) 


0.7952(66) 


32 


0.7562(91) 


0.7779(56) 


0.7812(220) 


0.7833(70) 


L — > oo 


0.720(16) 


0.764(14) 


0.735(28) 


0.744(17) 


XVd-o.f 


1.149/2 


0.208/1 


1.565/3 


0.025/1 


prob 


0.563 


0.649 


0.667 


0.874 



TABLE XI: Thermal exponents for the diluted 0(3) model 
with p=0.5. In the first column we have only used data with 
L > 8 while in the second and fourth ones we have only use 
data with L > 12. 
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0.055 
0.05 

0.045 
0.04 

0.035 
0.03 



: 1.0 

:0.9 
:0.7 
:0.5 







0.05 



0.1 



0.15 



0.2 



17 



FIG. 8: Joined extrapolation to L — *• oo for the r\ exponent de- 
duced from the vectorial susceptibility (x )■ Extrapolations 
are carried out by choosing a common value for the first term 
of Eq. (|22p for all dilutions, and by minimizing the joined \ 2 - 



0.05 



0.045 



0.04 



0.035 



0.03 



: 1.0 

:0.9 

= 0.7 
= 0.5 



o 



0.05 



0.1 



0.15 



17 



FIG. 9: Joined extrapolation with all p values to L — > oo 
for the rj exponent deduced from the vectorial magnetization 
{M v ). 




FIG. 10: Joined extrapolation with all p values to L 
for the v exponent deduced from the dpgX . 



0.95 



0.85 - 




0.75 - 



0.65 



FIG. 11: Joined extrapolation to L — > oo for the v exponent 
deduced from the dpg J. 



0.2 





V 


rfr 


X 


M 


XT 


M T 


Our results 


0.0390(9) 


0.0389(10) 


1.4251(13) 


1.4251(14) 


X 2 /d.o.f 


6.675/12 


5.104/15 


9.151/10 


13.931/11 


prob 


0.878 


0.991 


0.518 


0.237 


Ref. [11] 


0.0378(6) 









TABLE XII: Joined extrapolation with all p values for the 
magnetic exponent rj compared with the results from Ref. [ll| ] . 
The first three rows correspond to our L — > oo extrapolation, 
being prob the probability to find a larger value for the \ 2 °f 
the fit (it is a measure of the goodness of the fit) and d.o.f . 
the number of degrees of freedom (the total number of data 
minus the total number of adjustable parameters in the fitting 
function). 





V 






dfjgl 




Our results 


0.7126(46) 


0.7129(31) 


0.7294(81) 


0.7089(32) 


X 2 /d.o.f 


4.831/11 


6.606/6 


9.009/13 


9.609/7 


prob 


0.939 


0.359 


0.772 


0.212 


Ref. [11] 


0.7113(11) 









TABLE XIII: Joined extrapolation with all p values for the 
thermal exponent v compared with the results from Ref. [llj. 





gX 


T 


Our results 


0.62018(6) 


0.51366(19) 


X 2 /d.o.f 


10.324/9 


5.980/10 


prob 


0.325 


0.817 


Ref. [11] 


0.6202(1) 




Ref. [3] 


0.625783 





TABLE XIV: Joined extrapolation to L — * oo with all p values 
for the Binder cumulant, ,94, defined in Eq. (|lip . compared 
with results from Refs.fTJ and 0]. 
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